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Abstract 

In this paper, we first construct a model for transient free surface flows that takes into account the 
air entrainment by a sytem of 4 partial differential equations. We derive it by taking averaged values 
of gas and fluid velocities on the cross surface flow in the Euler equations (incompressible for the fluid 
and compressible for the gas). Then, we propose a mathematical kinetic interpretation of this system to 
finally construct a well-balanced kinetic scheme having the properties of conserving the still water steady 
state and possesing an energy. Finally, numerical tests on closed uniforms water pipes are performed and 
discussed. 
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■ 1 Introduction 

m 

The hydraulic transients of flow in pipes have been investigated by Hamam and McCorquodale [H] , Song et 
al. [19l [20] [21] , Wylie and Streeter [25] and others. However, the air flow in the pipeline and associated air 
pressure surge was not studied extensively. 

Two component gas-liquid mixture flows occur in piping systems in several industrial areas such as nuclear 
power plants, petroleum industries, geothermal power plants, pumping stations and sewage pipelines. Gas 
may be entrained in other liquid-carrying pipelines due to cavitation, [24) or gas release from solution due 
to a drop or increase in pressure. 

Unlike in a pure liquid in which the pressure wave velocity is constant, the wave velocity in a gas-liquid 
mixture varies with the pressure. Thus the main coefficients in the conservation equation of the momentum 
are pressure dependent and consequently the analysis of transients in the two-component flows is more 
complex. 

The most widely used analytical models for two-component fluid transients are the homogeneous model [3], 
the drift flux model [TOj [13l [14], and the separated flow model, [22l [23]. 

In the homogeneous model, the two phases are treated as a single pseudo-fluid with averaged properties 
[9] [25]: there is no relative motion or slip between the two phases. The governing equations are the same as 
the one for a single phase flow. However the inertial and gravitationnal effects can play a more important 
role than the relative velocity between the air phase and the liquid one: the interaction between the two 
phases has to be taken into account. 



On 
O 



- 1—1 

X 



*e-mail: Christian.Bourdarias@univ-savoie.fr 
t email: Mehmet.Ersoy@univ-savoie.fr, 
t email: Stephane.Gerbi@univ-savoie.fr 



1 



In the drift-flux model, [TU1 EH [H] , the velocity fields are expressed in terms of the mixture center-of-mass 
velocity and the drift velocity of the vapor phase, which is the vapor velocity with respect to the volume 
center of the mixture. The effects of thermal non-equilibrium are accommodated in the drift-flux model by 
a constitutive equation for phase change that specifies the rate of mass transfer per unit volume. Since the 
rates of mass and momentum transfer at the interfaces depend on the structure of two-phase flows, these 
constitutive equations for the drift velocity and the vapor generation are functions of flow regimes. 
Two-phase flows always involve some relative motion of one phase with respect to the other; therefore, a 
two-phase flow problem should be formulated in terms of two velocity fields. A general transient two- phase 
flow problem can be formulated by using a two-fluid model. The separated flow model considers the two 
phases individually, interacting with each other. Generally, this model will be written as a six equations 
partial differential system that represent the balance of mass, momentum and energy [221 123) . 
In the present paper, we consider the air entrainment appearing in the transient flow in closed pipes not 
completely filled: the liquid flow (as well as the air flow) is free surface. We thus derived the present model 
in a very closed manner as the derivation of Saint- Venant equations for free surface open channel flow: we 
consider that the liquid is incompressible and its pressure is only due to the gravitationnal effect. It is the 
well known hydrostatic pressure law. The air phase is supposed to be compressible and isothermal: thus 
we can use an equation of state of the form p = kp 1 where 7 represent the adiabatic exponent. To connect 
the two phases, we write the continuity of the normal stress tensor at the free surface separating the two 
phases: the hydrostatic pressure law for the fluid is thus coupled to the pressure law governing the air flow at 
the free surface. Then by taking averaged values over cross-sections of the main flow axis, we obtain a four 
equations partial differential system: we called it the two-layer water/air model. We state in theorem 12.11 
its main properties: the existence of a total energy despite the fact that the obtained system is not strictly 
hyperbolic (except for the particular case where the averaged speeds of the liquid and the gas are the same). 
The derivation of the proposed model and the investigation of the mathematical properties are presented in 
section [2] 

As we will see, the main difficulty in the construction of a numerical scheme to solve the obtained model 
comes from the fact that the proposed model is not hyperbolic: the eigenvalues of the jacobian of the flux 
may be complex. Moreover we cannot compute exactly these eigenvalues since we do not suppose that the 
velocity of each phase are close. This is the main reason why we chose to interpret this model as a kinetic 
model: we present the mathematical kinetic formulation of the four equations partial differential equations 
in the section [3J This mathematical formulation permits us to construct a well-balanced kinetic scheme 
presented in the section 0] 

Finally, we propose in section [5] some numerical simulation to focus on the role of the air entrainment on the 
free surface fluid flow. 

2 A two-layer water/air model 

We will consider throughout this paper a transient flow in a closed pipe that is not completely filled of water 
nor of air. The flow is composed by an air and water layer separated by a moving interface. We assume 
that the two layers are immiscibles. We will consider that the fluid is perfect and incompressible whereas 
the air is a perfect compressible gas: we neglect here the air-water interactions as condensation/evaporation 
and we focus only on the compressibility effects of the air on the water layer. The air layer is supposed to 
be isothermal and isentropic. The derivation of the two-layer model is then performed with the use of the 
3D Euler (incompressible for the fluid, compressible for the gas) equations. 

Remark 2.1 Throughout the paper, we use the following notations: we will denote by the index w the water 
layer and by the index a the air layer. We also use to the index a to write indifferently the liquid or gas 
layer. 

The pipe is assumed to have a symetry axis C which is represented by a straight line with, for the sake of 
simplicity, a constant angle 9 with the horizontal direction. (0,i, ,j,k) is a reference frame atached to this 
axis with k orthogonal to i in the vertical plane containg C, see Fig. [TJ We denote by VL t: a and flt,w the gas 
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and fluid domain occupied at time t. Then, at each point uj(x), we define the water section Q w (t, x) by the 
set: 

Q w (t,x) = {{y,z)eR 2 ;ze [—R(x),—R(x) + H w (t,x)], y G [0i(x,z),0 r (x,z)]} 
and the air one by the set: 

n a {t,x) = {(y,z) G K 2 ;ze [-R(x) + H w {t,x), R(x)], y G [/3 ; (x,z),/3 r (a;,z)]} 

where R(x) denotes the radius, H w (t,x) the water height at section Q w {t, x), H a (t,x) the air height at 
section Q a (t, x). I3i{x, z), f3 r {x, z) are respectively the left and right y values of the boundary points of the 
domain at altitude -R(x) < z < R(x) (see Fig. [Hand Fig. O 




Figure 1: Geometric characteristics of the domain 




Figure 2: Cross-section of the domain 



We have then the first natural coupling: 

H w (t,x)+H a {t,x) = 2R{x). (1) 
In what follows, we will denote h w (t,x) = —R(x) + H w (t,x) the algebraic water height. 
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2.1 The fluid free surface model 



The shallow water like equations for free surface flows are obtained from the incompressible Euler equations 
@ which writes in cartesian coordinates (x,y,z): 

div(U w ) =0, on R x Q, t , w , , 

9 t (U w ) + div(U w (g)U w ) + VP K , = F, onlxn,,, U 

where U w (t, x, y, z) = (U W ,V W ,W W ), P w (t,x,y,z), F denote respectively the fluid velocity, the scalar pressure 
and the exterior strength which represents the gravity. 

A kinematic law at the free surface and a no- leak condition U w .n w b = on the wetted boundary (n w b 
beeing the outward unit normal vector Fig. O complete the system (see the works of Bourdarias et al. for 
a complete derivation [5J [HI [7] ) . 

To model the effects of the air on the water, we write the continuity of the normal stress at the free boundary 
layer separated the two phases. Thus the water pressure P w may be writen as 

g(h w - z) cos 9 + P a (p)/p 

where the term g(h w — z) cos 9 represents the hydrostatic pressure whereas P a (jj) represents the air pressure, 
pit, x, y, z) — po represents the density of water at atmospheric pressure and ~p is the mean density of the air 
over the section O a (i, x). This is the second natural coupling between the two phases. 

Let us now introduce the conservative variables A(t, x) and Q(t, x) — A(t, x)u(t, x) representing respectively 
the wetted area and the discharge defined as: 



A(t,x) = / dydz, Q(t,x) — A(t,x)u(t,x) 
where u is the mean value of the speed 

u(t,x) = —^ — - f U w (t,x,y, z) dydz. 



A(t,x) 



We take averaged values along the cross-section Q w (t,x) in the incompressible Euler equations, and we 
approximate the averaged value of a product by the product of the averaged values ( see for instance 
[HISJIEIII])- This leads to the free surface model: 

f d t A + d x Q = 

dtQ + d x (%+ AP a (p)/p Q + gh(x, A) cos 9) = -gAd x Z 

\ A ) (3) 

+gl2(x,A) cos 9 

+P a (p)/ Po d x A 

where g is the gravity constant an Z(x) is the elevation of the point u)(x). The terms Ii(x,A) and l2(x,A) 
are defined by: 

Ii(x,A)= / (h w — Z)a(x, z) dz 

J-R 

and 

r h,„ 



I^ix, A) — / [h w — Z)d x a(x, z) dz . 

J-R 

They represent respectively the term of hydrostatic pressure and the pressure source term. In theses formulas 
cr(x, z) is the width of the cross-section at position x and at height z (see Fig. [2|). 
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2.2 The air layer model 

The derivation of the air layer model is based on the Euler compressible equations (0J which writes in 
cartesian coordinates (x,y, z): 

d t p a + div(p a U a ) = 0, Oil K X fl t ,a /,s 

d t (p a U a ) + div(p a U a ® U a ) + VP a = 0, on K x Q Ul [ > 

where U a (i, x, y, z) = (U a , V a , W a ) and p a (t, x, y, z) denotes the velocity and the density of the air whereas 
P a (t, x,y, z) is the scalar pressure. Let us mention here that we neglected the gravitationnal effect on the 
air layer. We neglect the second and third equation in the momentum equations to assume that the mean 
air flow follows the ir-axis. 

Remark 2.2 We can also obtain it from a classical rescaling e = 2R/L where R is the radius of the 
pipe section, L its length and we get in particular the equations d y p a = d z p a = 0. Thus in the following 
development, we may assume that the density of air depends only on t and x, see J2|/ for a formal asymptotic 
derivation of the model for mixed flows in closed pipes. 

We use the standard notation for the mean value of the density of air ~p over the air section fi a . We assume 
that the air layer is isentropic, isothermal, immiscible and follows a perfect gaz law. Thus we have the 
following equation of state: 

P a (p) = kp-< with A: = ^ (5) 

Pa 

for some reference pressure p a and density p a . The adiabatic index 7 is set to 7/5. 

With the definition of the free surface pressure at Subsection 12.11 and the air pressure definition above 
and assuming that the velocity is normal at the free surface of each layer, the two-layer model follows the 
continuity property of the normal stress for two immiscible perfect fluids at the interface. 
Let us now introduce the air area A by: 

A = dydz 

and the averaged air velocity v (we recall that u is the averaged velocity of the fluid layer) by: 



v ( t > x ) = ~jiu — \l U a {t,x,y,z)dydz . 
Ait, x) J na 

The conservative variables are M = ~p/poA and D = Mv and represent the rescaled air mass and discharge. 
Let us remark that we use rescaled mass and discharge instead of the real mass and discharge to be consistent 
with the free surface model ((3]). As done before to obtain the water layer model, we take averaged value 
in the Euler equations ((4J over sections Q a , and perform the same approximations on averaged values of a 
product to obtain the following model: 

d t M + d x D = / p/po (<9 t m + vd x m — v) .n ds 

d t D + d x (j^ + P a (p)/p A) = P a (p)/pod x (A) (6) 

+ / p/Po v {9 t m + vd x m — v) .nds 

where v is the velocity in the (j, k)-plane. The boundary dfl a is divided into the free surface boundary T f s 

m 

and the border of the pipe in contact with air r c . For m S T c , n = - — - is the outward unit vector at 

H 

the point m in the f2-plane and m stands for the vector umn while n denotes -nj s on Tj s . As the pipe is 
infinitely rigid, the non-penetration condition holds on r c , namely: U a .n a ft = (see Fig. [1]). 
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On the free surface, as the boundary kinematic condition used for the water layer at the free surface holds, 
the integral appearing in System ([6]) vanishes. Using the equation of state ([5]), and using the air sound speed 
defined by 



we get finaly the air-layer model: 



7-1 



d t D + d x 



dp 



V A J 



d t M + d x D 

D 2 M 





M 

1 



cld x (A) 



(7) 



(8) 



2.3 The two-layer model 

The two-layer model is then simply obtain by the apposition of the models © and ([5]) and using the coupling 
([T]) which also writes A + A = S where S — S(x) denotes the pipe section: 



dtD + d x 



8 t M + d x D 
D 2 M 
~M + ~' 





M 
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c 2 a d x {S-A) 



d t A + d x Q = 



d t Q + d x [ ^-+gI 1 {x,A)cos9 + A 



(9) 



7(S-A) 



~gAd x Z + gl 2 (x, A) cos 8 
clM 



l(S-A) 



d x A 



The two-layer model has the inconvenient to be non hyperbolic since we may have complex eigenvalues 
as we will see. Moreover, we have no way to give it explicitly. Indeed, setting the conservative variable 
W = (M, D, A, QY, let us consider System ([9]) with 8 = and S(x) = S (for the sake of simplicity) written 
under the non conservative form: 

(10) 



9 t W + V(W)d x W = TS(W) 



with the convection matrix: 



V 



A 



\ (S-A) a 



1 

2v 







M 2 

S~~A° a 


AM 2 
{S-A) 2 ° a 



\ 

1 

2u j 



where c a being the air sound speed defined by ([7|) and c„ 
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T(A) 



cos 8 the water sound speed. T(A) is 



the width of the free surface at height h w (A). 

Since the angle 8 is set to and S is constant, the source term TS(W) vanishes. The eigenvalues A of the 
convection matrix must satisfy the following polynomial equation: 



(A 2 - 2v\ + (c 2 a - v 2 )) (V - 2u\ - [c 2 w + 



AM 



AM 
JS^AJ 2 



(i 



Remark 2.3 We can remark that when the speed of each layer are equal or even close, we can write, explicitly 
in the case u — v or an asymptotic development in the case where u w v, the four distinct eigenvalues: 



1 / , / AM \ „ \ f n f AM \ „ V AM 

Two of them are real and two are complex. This makes the system non-hyperbolic. Thus the classical Roe like 
Finite Volume method is not appropriate since it is necessary to solve a Riemann problem at the interfaces 
of the cells which asks to compute the eigenvalues. 

Moreover, if we assume that the air density is small enough (e.g. c a « 0), the eigenvalues are simply 
v, v, u± Cyj. However, in the case of a high speed hydraulic jump, the air speed can reach severe values 
and may be very different of the speed of the water. Also the air layer may be strongly compressed and the 
eigenvalues cannot be computed. 

Nonetheless, the problem of complex eigenvalues can be solved easily in a kinetic framework (see sections^ 
and^4ty. 

Although System ([9|) is not hyperbolic when u is different from v, we may construct a non increasing total 
energy, which permits us to think that there exist entropic solutions of System |9]) (see also [3]). The results 
are summarized in the following theorem: 

Theorem 2.1 

1. For smooth solutions, the velocity u and v satisfies 

^ + ^(y + ^tt) = ( n ) 

(u 2 c 2 M \ 

d t u + d x ( — + gh w (A) cos 9 + gZ + ° _ j = . (12) 

2. System (0) admits a mathematical total energy 

E = E a + E w 

with 

Mv 2 c 2 M 

&a, = — 7 



7(7 - 1) 

and 

Au 2 

E w = — + gA(h w - I x (x, A) j A) cos 6 + gAZ) 
which satisfies the entropy equality 

d t E + d x H = 
where the flux H is the total entropy flux defined by 

H = H a + H w 

where the air entropy flux reads: 



and the water entropy flux reads 

H w = [E w + gh (x, A) cos9 + A 



cjM 
(S-A) 
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3. The energy E a and E w satisfies the following entropy flux equalities: 

d t E a + d x H a = . C ° — d t A (13) 
7{S-A) 

and 

c 2 M 

d t E w + 8 x H w = j§ -dtA . (14) 

7(6 -A) u 

Proof of theorem 12.11 

The proof of these results relies only on algebraic computations and combinations of the equations forming 
System ©. 

Remark 2.4 Let us remark that Equations (13\) - [L{\ ) can be interpreted in such a way: the system is physi- 
cally closed e.g. the energy dissipated (or gained) by the water layer is forwarded to (or lost by) the air layer. 
This property ensures the total energy equality dtE + d x H = 0. 

3 The mathematical kinetic formulation 

As pointed out before (see Subsection 12. 3p . Roe like Finite Volume method is not appropriate since, in our 
case, the eigenvalues cannot be computed explicitly and approximations cannot be used unlike the multilayer 
Saint- Venant system solved by a Roe like Finite Volume method in [2]. 

A kinetic approach is proposed to overcome this difficulty We recall the kinetic formulation (see e.g. [T7] ) 
and we apply it to the two-layer model ([9]) . Let us consider a smooth real function \ which has the following 
properties: 

= > , / xHdio = 1, f Lu 2 x(Lu)aio = 1 . (15) 

We define the Gibbs equilibrium as 

b a \ b a 

which represents the density of particles at time t, position x with the kinetic speed £ with 

r2 



M 



h(x,A) c\M 
g __ C0S fl + ___ lf 



We denote: 



. if a — w n J Q if a — w , J u if a — w 

a = { M if a = a ' Qa = \ D if a = a Ua ~ \ v if a = a ' 

Thus we may state the following kinetic interpretation of System © : 

Theorem 3.1 The couple of functions (A a , Q a ) is a strong solution of System (0) if and only if M. a satisfies 
the two kinetic transport equations (one per phase): 

dt-M a + ^dxM. a + (ftad^Ma = K a (t, x, £) for a — a and a = w (16) 

for some collision term K a (t,x,£) which satisfies for a.e. (t,x) 

K a d£ = 0, [tK a dt = . (17) 
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The function <p a is defined by: 



^-d x {S-A) if a = A 



S — A {1-8) 
^Z-g^^-cose-^jcldMA) if a = W 

Proof of theorem 13.11 The proof relies on a very obvious computation. Indeed, the equations forming the 
System (J9)) are obtained by taking the moments of the two kinetic Equations (fTBj) with respect to £ against 
1, £ and £ 2 : the right-hand side vanishes according to (JTTJ) and the left-hand sides coincides exactly thanks 
to (|15p . These are consequences of the following relations verified by the microscopic equilibrium: 



A a = M„(M,^ 
Js, 

Q a = J tM a {t,x,£)d£ 

+ b 2 a A a = f ?M a (t,x,S)dS 
Jn 



Remark 3.1 The previous results \3.1\) are obtained assuming a constant angle 9. The generalization of 
these results with varying = 6{x) can be found in JS$ for the free surface layer. 



4 The kinetic scheme 

This section is devoted to the construction of the numerical kinetic scheme and its properties. The numerical 
scheme is obtained by using the previous kinetic formulation (|16p . The term <f> a is easily upwinded. Indeed, 
the conservative term d x Z is a given piecewise constant function and all n on conservative products of 4> a 
are defined following, for instance [11[|15[[16] . For the sake of simplicity, we will only consider uniform pipes, 
that induces S{x) = S. 

In the sequel, for the sake of simplicity, we consider infinite space domain. Let us consider the discretisation 
{rrii)i£z of the spatial domain with 

mi = {x^i/2, x i+1 / 2 ), hi — x i+ i/2 — ^i-i/2) « G Z 

which are respectively the cell and mesh size. Let At" = t n+ i — t n , n <E N be the timestep. 

Q n ■ 

Let W"j = {A™ i,Qa i), u™ i = — ^- be respectively the approximation of the mean value of {A a , Q a ) and 

A n ( £ — u n - \ 

the velocity u a on m% at time t„. Let A-l™ i (£) = \ — , n "'^ De the approximation of the microscopic 

a, i \ a,i J 

quantities and </>™ jl mi (X) of <f> a . Then, integrating System © over ]xi_i/ 2 ,x i+ i[x[t n ,t n+ i[, we can write 
a Finite Volume scheme as follows: 

W" +1 = W" ■ - — - I F , - F + , 1 (19) 

vv q,i vv q,i ^ \ r a,i+l/2 r a,i-l/2J- \ la ) 

where the numerical flux F^ i±1 ^ 2 are defined from the microscopic fluxes M^ i±1 / 2 with the relation: 
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The microscopic fluxes Mai±i/2 are upwinded by the following formulas as done in [¥1 IS1 ITS]: 

positive transmission reflection 



M-, i+1/2 (t) = l {( >o } M n aA (0 +t {i<042 _ 2gA<f , ai+y2<0} Ml i (-0 

+ l{ € <o,^_ 2fl A^ , i+1/2 >o} M a,i+i (- \/e 2 -2.gA0 Q:l+1/2 ) 
" , ' 

negative transmission 
negative transmission reflection 



(20) 



*£,i+i/a(0 = 1 {f<o}^2, i+ i(0 + 1 { e >o, f = + 2 fl A^, <+I/2 <o}^2, i +i(-0 
+ 1 {^>o,^+2 S A^, i+1/2 >o}>'S,i ( v /e 2 + 2ffA^ Qil+1/2 ^ 



with 



A<f>a,i+l/2 = < 



positive transmission 



^ (A i+1 - ^) if 

7 



r 2 M \ 1 

.</ ( Z; , - Z, ) - ( - HA Z+1 /A t ) if a u; 



To construct this upwinding fluxes, we used firstly the fact that every quantities involved are piecewise 
constant functions and may define non conservative products. Using the "straight lines" paths: 

$(s, Wi, W l+1 ) = sW i+ i + (1 - s)Wi, s e [0, 1] 

(see e.g. [Tl]|T5j[T6] and Fig. [3j with Wi, Wi+i the left and right state at the discontinuity Xi+i/% permits 
us to approach any non conservative product f(x, W)d x W as /(Wj+i — Wi) with the notation 

/= / f(x i+ i/2,$(s,Wi,W i+1 ))ds . 
Jo 

A^ a ,i±i/2 denotes the upwinded source terms corresponding to (|18p and £ 2 ±2A0 Q 4+1/2 is the jump condition 
for a particle with the kinetic speed £ which is necessary to: 

• be reflected: this means that the particle has not enough kinetic energy £, 2 /2 to overpass the potential 
barrer which corresponds to potential energy A(f> ai 

• overpass the potential barer with a positive speed, 

• overpass the potential barer with a negative speed, 
as displayed on Fig. [31 
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Reflection & Transmission 

t> 

<t ' <t>i+i 




Reflection 

Positive transmission 
Negative transmission 
Reflection zone 



M?({) X ' +1 ' 2 -M? +I R) 

Figure 3: The potential barer: transmission and reflection of particle 
Top: the physical configuration 
Middle: the characteristic solution in (X, H)-plane 
Bottom: the characteristic solution in (X, t)-plane 

We use as \ function the function defined by xi w ) — — 7=^1-^ J3]( w ) ■ Thus a stability CFL condition is 

2V3 ' 

necessary and under this condition, we can state the following results on the positivity of the air and water 
areas, the still steady states and the flooding and drying process: 



Theorem 4.1 Assuming the CFL condition 



hi 



At n < min . 

iez Imax^K^I + v^) 



holds, we have: 

1. the numerical kinetic scheme keeps A a positive. 

2. The still air/water steady state is preserved: 

u? = , v? = 
Zi 



„n 2 



MP 



9l(Sr - A?) 



est i 



cst 2 



3. Drying and flooding areas are also treated. 
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5 Numerical tests 



The numerical validation consist in showing the air entrainment effect in closed pipes with constant section 
and no slope. Since experimental data for such flows in any pipes are not available, we focus only the 
qualitative behavior of our model and our method for different upsrcam and downstream condition. 
The numerical experiments are performed in the case of a horizontal 100 m long closed uniform circular pipe 
with diameter 2 to. 

The simulation starts from the inital steady state which consists of a free surface flow with a null discharge 
Q w = m 3 /s and a water height of y w = 0.2 to. 

In the following tests, the upstream water height is increasing in 25 s at the height y w — 0.4 m and the air 
density at the downstream boundary is taken successively as p a = 1.29349 kg/m 3 and 1.29349 10 -2 kg/m 3 . 
The downstream discharge of water and air is kept equal 0to 3 /s (as displayed on Fig. 13} . 
We compare then the results obtained for the present model and the Saint- Venant Equations for the water 
layer only, solved by the kinetic scheme at the middle of the pipe for the water heigh, the water discharge, 
the air pressure. The piezometric head is defined as: 

piezo = z + h w with h w the water height 




(a) Downstream Air and Water discharge (b) Upstream Water piezomeric head (c) Upstream Air discharge 

Figure 4: Upstream and downstream boundary conditions 



i ... 




(a) Single fluid (b) Piezometric head of water (c) Air pressure 

Figure 5: Piezometric head of water and air pressure at the middle of the pipe. p a = 1.29349 10 -2 kg/m 3 
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(a) Single fluid (b) Water discharge (c) Air discharge 

Figure 6: Discharge of water and air at the middle of the pipe. p a — 1.29349 1CP 2 kg/m 3 





(a) Single fluid (b) Piezometric head of water (c) Air pressure 

Figure 7: Piezometric head of water and air pressure at the middle of the pipe. p a — 1.29349 kg/m 3 




(a) Single fluid 



(b) Water discharge 



(c) Air discharge 



Figure 8: Discharge of water and air at the middle of the pipe. p a = 1.29349 kg/m 3 

We can see on Fig. [5] and Fig. [6] (at the downstream and upstream boundaries p a = 1.29349 10 -2 kg/m 3 ) 
that the influence on the air entrainment on the water flow is negligeable whereas on Fig. [7] and Fig. [8] 
(at the downstream and upstream boundaries p a — 1.29349 kg/m 3 ) the air entrainment is very important 
: since we imposed that the discharge of air is null at the two boundaries, the air layer is pushed in the 
pipe, so that the air layer is compressed. On the other hand, the progression of discharge of the water is 
accelerated since the air layer pushes on the water layer. 

6 Conclusion 

This study is of course the first step in the comprehension and the modelisation of the air role in the transient 
flows in closed pipes. We have now to treat the air entrapment pocket which may be encounter in a rapid 
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filling process of a closed pipe: in this case a portion of the pipe will be completely filled and air pocket 
will be entrapped. This is the next step of our research work since in previous works, we have derived a 
model for pressurised flow or mixed flows in closed pipes without taking into account the role of the air and 
proposed a Roe like Finite Volume method and a kinetic scheme [H [51 G3 El [5] . 

Next, we have to deal with the evaporation/condesation of gas or water when it is not supposed to be 
isothermal and lately, because in the water hammer phenomenon, large depression may occur, we have to 
deal with the natural cavitation problem. 
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